import numpy as np  
import matplotlib.pyplot as plt  
from mpl_toolkits.mplot3d import Axes3D  
  
# 洛伦兹系统的参数  
sigma = 10.0  
rho = 28.0  
beta = 8.0 / 3.0  
  
# 时间步长和总时间  
dt = 0.01  
total_time = 50  
  
# 初始条件  
x0, y0, z0 = 0.0, 1.0, 1.05  
  
# 时间数组  
t = np.arange(0, total_time, dt)  
  
# 初始化变量数组  
x = np.zeros(len(t))  
y = np.zeros(len(t))  
z = np.zeros(len(t))  
  
# 设置初始值  
x[0], y[0], z[0] = x0, y0, z0  
  
# 迭代计算洛伦兹系统  
for i in range(1, len(t)):  
    dx = sigma * (y[i-1] - x[i-1])  
    dy = x[i-1] * (rho - z[i-1]) - y[i-1]  
    dz = x[i-1] * y[i-1] - beta * z[i-1]  
      
    x[i] = x[i-1] + dx * dt  
    y[i] = y[i-1] + dy * dt  
    z[i] = z[i-1] + dz * dt  
  
# 绘制洛伦兹吸引子  
fig = plt.figure()  
ax = fig.add_subplot(111, projection='3d')  
ax.plot(x, y, z, lw=0.5)  
ax.set_xlabel("X Axis")  
ax.set_ylabel("Y Axis")  
ax.set_zlabel("Z Axis")  
ax.set_title("Lorenz Attractor")  
plt.show()